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ABSTRACT 

The influence of two-dimensional finite element modeling assumptions on the debonding prediction 
for skin-stiffener specimens was investigated. Geometrically nonlinear finite element analyses using 
two-dimensional plane-stress and plane-strain elements as well as three different generalized plane strain 
type approaches were performed. The computed skin and flange strains, transverse tensile stresses and 
energy release rates were compared to results obtained from three-dimensional simulations. The study 
showed that for strains and energy release rate computations the generalized plane strain assumptions 
yielded results closest to the full three-dimensional analysis. For computed transverse tensile stresses the 
plane stress assumption gave the best agreement. Based on this study it is recommended that results from 
plane stress and plane strain models be used as upper and lower bounds. The results from generalized 
plane strain models fall between the results obtained from plane stress and plane strain models. Two- 
dimensional models may also be used to qualitatively evaluate the stress distribution in a ply and the 
variation of energy release rates and mixed mode ratios with delamination length. For more accurate 
predictions, however, a three-dimensional analysis is required. 


BACKGROUND 

Many composite components in aerospace structures are made of flat or curved panels with co- 
cured or adhesively bonded frames and stiffeners. Previous investigations of the failure of secondary 
bonded structures focused on loading conditions typically experienced by aircraft crown fuselage panels. 
Tests were conducted with specimens cut from a full-size panel to verify the integrity of the bondline 
between the skin and the flange or frame [1], A simpler and cheaper specimen configuration that would 
allow detailed observations of the failure mechanism at the skin/flange interface was proposed in reference 
[2]. The investigations focused on the failure mechanisms of a bonded skin/flange coupon configuration 
subjected to tension, three and four-point bending, and combined tension/bending loading [3-6]. 
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An analytical methodology was also developed to predict the location and orientation of the first 
transverse matrix crack based on the principal transverse tension stress distribution in the off axis plies 
nearest the bondline in the vicinity of the flange tip [7]. Earlier investigations [2, 3] indicated that the 
matrix cracking occurred when the maximum principal transverse tensile stress, G tt , normal to the fiber 
direction, as calculated by 
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reached the transverse tension strength of the material. Here G 22 and G 33 denote the normal stresses and 
T 23 the shear stress in the 2-3 plane perpendicular to the fiber direction. 

A fracture mechanics approach was used to investigate delamination onset once the initial crack had 
formed. The initial crack was modeled as a discrete discontinuity at a location suggested by the 
microscopic investigation. The Virtual Crack Closure Technique (VCCT) was used to calculate mixed 
mode strain energy release rates G l (due to crack opening), G n , (due to in-plane shear) and G m (due to in- 
plane scissoring) [8-10]. Computed total strain energy release rates, G T , were compared to the critical 
value, G c , of the material for the appropriate mixed mode ratio (G n /G T ) in order to determine the potential 
for delamination growth [7]. 


INTRODUCTION 

The current study focused on skin-stiffener debonding resulting from buckling of a thin-gage 
composite fuselage structure as described in reference [7]. The specimens consisted of a bonded skin and 
flange assembly as shown in Figure 1(a). An IM7/8552 graphite/epoxy system was used for both the skin 
and flange. The skin was made of prepreg tape with a measured average ply thickness of h =0.148 mm 
and had a [45/-45/0/-45/45/90/90/-45/45/0/45/-45] lay-up. The flange was made of a plain-weave fabric 
with a thickness of h =0.212 mm. The flange lay-up was [45/0/45/0/45/0/45/0/45]^ where the subscript 
“f ’ denotes fabric, “0” represents a 0°-90° fabric ply and “45” represents a 0°-90° fabric ply rotated by 
45°. The measured bondline thickness averaged 0.178mm. Specimens were 25.4-mm wide and 177.8-mm 
long. The properties of the graphite/epoxy material and the adhesive were measured at Boeing and are part 
of the standard design database for the V-22 tilt-rotor aircraft. Typical properties are summarized in 
Table 1. 

Four quasi-static tension tests with a gage length of 101.6 mm were performed as outlined in [7]. 
The value of the damage onset load was averaged from four tests and determined to be P =17.8 kN with a 
coefficient of variation of 8.9%. The tests were terminated when the flange debonded from the skin. 
Damage was documented from photographs of the polished specimen edges at each of the four flange 
corners identified in Figure 1(a). Comers 1 and 4 and corners 2 and 3 had identical damage patterns. The 
damage at corners 2 and 3, formed first and consisted of a matrix crack in the 45° skin surface ply and a 
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delamination at the +457-45° interface as shown in Figure 1(b). Therefore, this damage pattern was the 
focus of the current and related earlier analyses [7]. 

The complex nature of the failure observed during the experiments, where the delamination changed 
across the specimen width from a delamination running at the skin surface 457-45° layer interface to a 
delamination propagating in the bondline above, suggests the need for a three-dimensional model [7]. 
Since many layers of brick elements through the thickness would be required to model the individual plies, 
the size of three-dimensional finite element models, however, may become prohibitively large. Two- 
dimensional models of a longitudinal cut through the specimen allow for a detailed modeling of the 
individual plies and the adhesive in thickness direction. Two-dimensional models are preferred by industry 
due to the fact that modeling time, as well as computational time, remains affordable, especially if many 
different configurations have to be analyzed during the initial design phase. These models have been used 
extensively during previous studies [2-4, 6, 7]. However, it is inherent to any two-dimensional finite 
element model that the geometry, boundary conditions and other properties are constant across the entire 
width. 

The current study presents an intermediate step where three-dimensional models were created by 
extruding two-dimensional models across the width. The fact that the delamination changed across the 
specimen width from a delamination running at the skin surface 457-45° layer interface to a delamination 
propagating in the bondline above, however, is still not accounted for in this model. Nevertheless, the 
three-dimensional model takes width effects into account and therefore provides insight into the limitations 
of the use of two-dimensional finite element models. For future detailed modeling and analysis of the 
damage observed during the experiments, the shell/3D modeling technique [11] offers great potential for 
saving modeling and computational effort because only a relatively small section in the vicinity of the 
delamination front needs to modeled with solid elements. 

The objective of the current study was to evaluate how the assumptions made in developing two- 
dimensional finite element models of the specimen would effect the damage onset prediction for skin- 
stiffener debonding. Geometrically nonlinear finite element analyses were performed using ABAQUS® 
two-dimensional plane-stress, plane-strain and generalized plane strain elements. Additionally, two 
different methods to create generalized plane strain conditions were studied [12, 13]. The results were 
compared to data obtained from full three-dimensional simulations to study the feasibility and limitations 
of using simple two-dimensional models. Computed transverse tensile stresses in the skin surface 45° 
layer, flange and skin strains, as well as mixed mode strain energy release rates were compared. 
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ANALYSIS FORMULATION 


Two-Dimensional Finite Element Models 

Several assumptions are possible in two-dimensional models. A plane-stress model imposes the out 
of plane stresses to be zero (G ZZ =T XZ = T yz =0) at the free edge of each ply and allows the displacement to 
be the free parameter. The plane-strain model, on the other hand, imposes the out of plane strains to be 
zero (£ zz =y x2 = Y yz =0) at the free edge, which excessively constrains the plies. Hence, the two dimensional 
plane stress and plane strain conditions may serve as upper and lower limits compared to the three- 
dimensional solution. 

An outline of the two-dimensional model of the specimen showing the boundary conditions, and 
loads applied is shown in Figure 2. Finite element models for an undamaged and a damaged specimen 
were developed using a refined mesh in the critical area of the 45° skin ply where cracking was observed 
during the tests. Outside the mesh refinement area, all plies were modeled with one element through the 
ply thickness. In the refined region, four elements were used per ply thickness for the first two individual 
flange plies above the bondline and the top skin ply below the bondline as shown in Figures 2(b) and (c) 
and described in more detail in reference [7] . Four elements through the thickness were also used to 
model the adhesive film. Ply properties and adhesive material properties used in this study are summarized 
in Table 1. 

Damage in the form of a matrix crack in the top 45° skin ply that developed into a delamination in 
the 457-45° interface is shown in Figure 1(b). Delaminations of various lengths were modeled by discrete 
discontinuities as shown in Figure 2(c). The delamination propagating in a bimaterial interface (457-45°) 
may, for small element lengths A a at the crack tip, result in incorrect mode separation of the total energy 
release rate. This effect is caused by the oscillatory singularity at the crack tip [14]. To avoid the effect, for 
the current investigation the element length A a was chosen to be about 1/3 of the ply thickness. At the 
opposite taper, the original mesh created for the model of the undamaged specimen was used. Finite 
element solutions were obtained using the commercial ABAQUS®/Standard finite element software. 
Eight-noded quadrilateral plane-stress (CPS8R) and plane-strain (CPE8R) elements with quadratic shape 
functions and a reduced (2x2) integration scheme were utilized for the geometric nonlinear analyses [15]. 

Generalized Plane Strain Models 

One alternative to a full three-dimensional simulation that requires about the same modeling effort 
as the simple two-dimensional models is a generalized plane-strain model. Generalized plane strain 
elements provided by ABAQUS® code are typically used to model a section of a long structure that is free 
to expand or is subjected to loading perpendicular to the plane of modeling. The formulation involves a 
model that lies between two planes that can move with respect to each other, and hence, cause strain in the 
direction perpendicular to the plane of the model that varies linearly with respect to position in the planes 
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[15]. Ten-noded quadrilateral generalized plane-strain (CGPE10R) elements with quadratic shape 
functions and a reduced (2x2) integration scheme were used in the current models. The models shown in 
Figure 2 were modified to use generalized plane strain analysis. Each element was assigned two additional 
nodes, which determined the position of the bounding planes and were common to all the elements. 

A different approach to create a generalized plane strain case was suggested by Minguet [12]. The 
generalized plane-strain case is regarded as an intermediate state between plane strain and plane stress, and 
assumes that the out of plane normal strain £ ZZ =-V L £ XX , where V L is the laminate Poisson’s ratio, and the 
shear strain y xz = Y yz =0. With these assumptions, ply stiffnesses were calculated for each ply angle in the 
specimen and entered as material properties for the two-dimensional models shown in Figure 3. 
Eight-noded quadrilateral plane-stress (CPS8R) elements with quadratic shape functions and a reduced 
(2x2) integration scheme were used for this analysis. 

Another possibility to create a two-dimensional generalized plane strain model was described by 
Konig [13]. The model is obtained by using one row of three-dimensional brick elements instead of plate, 
shell or membrane elements. Therefore, the two-dimensional models shown in Figures 2(b) and (c) were 
extruded into the models shown in Figures 3(b) and (c). This resulted in a 0.1 mm wide model made of 
one row of brick elements (C3D20R) with quadratic shape functions and a reduced integration scheme. 
The advantage of this model is the inclusion of three-dimensional effects in this two-dimensional finite 
element model. However, the number of degrees of freedom almost quadruple in comparison with 
membrane elements. The model may be compared to a wall, where one side of the model represents the 
midplane and the other side is a parallel plane. The midplane of the specimen as shown Figure 3(b) is a 
plane of symmetry and remains plane under the applied load. The displacements perpendicular to the 
plane are surpressed (vv=0) at all nodes in the midplane. All displacements in the z- direction at the nodes 
of the parallel plane are coupled to a constant value w. In ABAQUS® the coupling is achieved by grouping 
all nodes in the parallel plane except for one. The w displacement of all the nodes grouped is then set 
equal to the w displacement of the one node not included in the group using a linear constraint equation. 
These coupled displacements at the parallel plane will generate stress resultants at each node of this plane, 
but the summation of these stress resultants has to vanish because of the zero load in this direction. 

Three-Dimensional Finite Element Model 

The three-dimensional model of the specimen with load and boundary conditions is shown in 
Figure 4(a). The specimen was modeled with ABAQUS® solid twenty-noded hexahedral elements 
(C3D20R) with quadratic shape functions and a reduced integration scheme. Two-dimensional models of 
the undamaged and damaged specimen were created first in the x-y plane, and were then extruded into ten 
elements across the width of the specimen as shown in Figures 4(b) and (c). Although these three- 
dimensional models are not fine enough in the vicinity of the free edges ("=().() and z=25.4) to accurately 
capture the influence of the free edges on the distribution of the stresses and energy release rates across 
the width they are able to account for width effects that two-dimensional models do not include. 
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The same refined mesh as in the two-dimensional model (Figure 2) was used in the critical area of 
the 45° skin ply where cracking was observed during the tests. Outside the refined area, the mesh was 
modified to prevent the three-dimensional model from becoming excessively large (Figure 4(b)). The s ki n 
plies were grouped into three layered elements with -45/45/90, 90/-45/45 and 0/45/-45 respectively, thus 
taking advantage of the composite solid element option in ABAQUS® [15]. The fabric layers were 
grouped into two layered elements as shown in Figure 4(b). In the transition region at the flange several 
plies were modeled by one element with material properties smeared using the rule of mixtures [16]. The 
procedure did not ensure the full A-B-D contribution of the plies but appeared suitable for the small 
transition region to enforce a reasonable model size. The final model had 32,090 elements and 145,432 
nodes, yielding a total of 436,296 degrees of freedom. 


ANALYTICAL INVESTIGATION 
Global Response 

First, the global response of the specimens was computed at the mean quasi-static damage onset load 
(P=17.8 kN) determined from experiments. The load-strain behavior computed from different FE models 
was compared to the corresponding experimental results. This global response was used to examine 
whether the FE models, the boundary conditions, the loads and the material properties used in the model 
yielded reasonable results. Strains were averaged from computed nodal point values over a length 
corresponding to the dimensions of the strain gages shown in Figure 1(a) [7]. 

A schematic of the deformed geometry, the boundary conditions, and the load applied in the 
simulations is shown in Figure 2(a). In the schematic, the elongation of the specimen caused by the 
applied tensile load is shown. The bending effect caused by the load eccentricity in the flange region, the 
asymmetric layup with respect to the neutral axis, and the membrane stiffening effect is also shown. 

A comparison of measured strains at the surface of the flange and computed results is shown in 
Figure 5. The plane-strain model was stiffer, yielding an upper bound, and the plane-stress model was 
more compliant, yielding a lower bound compared to the results from the other models. The results from 
the generalized plain strain models and the full three-dimensional model agree best with the strain 
recorded during the tests. 

In Figure 6, measured strains at the surface of the top 45° skin ply near the flange tip (see Figure 1(a)) 
and computed surface strains were compared. Scatter in the experiments can not be shown as only one test 
was equipped with a strain gage on the skin. As mentioned above, the plane-strain models showed a stiffer 
behavior and the plane-stress models a more compliant behavior compared to the tested specimen. The 
results from two of the generalized plain strain cases and the full three-dimensional model fall between the 
results from the plane stress and plane strain simulation, however, appear stiffer when compared to this 
single test. The generalized plain strain model suggested by Konig [13]. yields results similar to the plane 
stress model. 
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Local Response 


A stress analysis was used to study the initial damage in the form of matrix cracking. The maximum 
principal tensile stress, O tt , normal to the fiber direction was calculated for the static failure load using the 
two-dimensional and three-dimensional finite element models introduced in the previous section. The 
stress distribution in the top 45° skin ply is plotted in the immediate vicinity of the flange tip in Figure 7. 
In the graph, a=0 corresponds to left grip as shown in Figures 2 to 4, and a- 2 6 . 4 corresponds to the left 
flange tip. During the fatigue tests, first matrix cracking was observed at locations right next to the flange 
tip which corresponds to the location where peak stresses were calculated [7]. 

The transverse tensile strength and scatter band, for IM7/8552, obtained from three -point bending 
tests of 90° lamina, were added to the plot for comparison [17]. Stresses at the static failure load computed 
from the plane-stress analysis correspond with the transverse tensile strength. Stresses from the plane- 
strain analysis are excessively high. This is caused by the constraints inherent to the plane-strain model, 
particularly in the ±45° plies. 

A fracture mechanics approach was used to investigate delamination onset once the initial crack had 
formed. During a series of nonlinear finite element analyses, strain energy release rates were computed at 
each front location for the loads applied in the experiments. A critical energy release rate, G c , needs to be 
determined to predict delamination onset. This critical G is generally identified based on the shape of the 
total energy release rate versus delamination length curve, which is determined through analysis as shown 
in Figure 8. The G T versus x curve reached a peak at some virtual delamination length and then decreased. 
The delamination was extended to a total simulated length of 2.2 mm to ascertain that the peak value had 
been captured. The total energy release rates computed for all two-dimensional and generalized plane 
strain models are plotted in Figure 8. The values obtained from three-dimensional analysis along the 
centerline of the specimen (z=12.7) was included in Figure 8. Qualitatively, all results follow the same 
trend. After a small initial drop the computed total energy release rate increases sharply with delamination 
length, reaches a peak value and gradually decreases. As expected the values from plane stress and plane 
strain analysis form an upper and lower bound, except for very short delamination lengths. Results from 
generalized plane strain models and the three-dimensional model fall between the results from two- 
dimensional analysis. The generalized plain strain model suggested by Konig yields results closer to those 
obtained from the plane stress model. All generalized plane strain peak G T results are within 8.5% of the 
results obtained from full three-dimensional analysis at the center of the specimen. 

The variation of mixed mode ratio G s /G, with delamination length is shown in Figure 9. Here G s 
denotes the sum of the in-plane shearing components G u + G np and G r denotes the total energy release rate 
Gj+G n +G ni , where G 7 is the opening mode. For two-dimensional analyses, where G n j= 0, this definition is 
equal to the previously used definition of the mixed mode ratio, G n /G r For three-dimensional analysis, 
which also yields results for the scissoring mode G np the modified definition of G s is introduced. All 
analyses yield the same trend, where the delamination initially starts with high shearing components, 
followed by a drop which is equivalent to an increase in G r For longer delaminations a gradual increase in 
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shearing components is observed. For a more realistic comparison of the G u contribution from two- 
dimensional and three-dimensional models results from the three-dimensional analysis were calculated 
assuming G s =G n and ignoring G m As Figure 9 shows, these results fall between the results from two- 
dimensional analysis, for longer delamination lengths. For G s =G n +G np the contribution of the in-plane 
shearing component is larger compared to the mode ratios obtained from two-dimensional analysis. 


CONCLUDING REMARKS 

The influence of the assumptions made in developing two-dimensional finite element models on 
skin-stiffener debonding specimens was studied. Geometrically nonlinear finite element analyses using 
two-dimensional plane-stress and plane-strain elements as well as three different generalized plane strain 
type approaches were performed. The computed skin and flange strains, transverse tensile stresses and 
energy release rates were compared to results obtained from three-dimensional simulations. 

The plane stress and plane strain models provided results for skin and flange strains, as well as 
energy release rates, which form an upper and lower bound of the results obtained from full three- 
dimensional analysis. The results from generalized plane strain models fall between the results obtained 
from plane stress and plane strain models. The generalized plane strain models capture the surface strains 
on the flange and skin very well and computed energy release rates are within 9% of the results from 
three-dimensional analysis. Transverse tensile stresses obtained from plane strain and generalized plane 
strain models, however were excessively high when compared to results from three-dimensional analysis. 
The stresses from full three-dimensional analysis were lower than results from any of the two-dimensional 
models, but are closest to the plane stress results. 

In general, two-dimensional models are preferred by industry because modeling time, as well as 
computational time, remains affordable, especially if many different configurations have to be analyzed 
during the initial design phase. Based on the results of this investigation it is therefore recommended to 
use results from plane stress and plane strain models as upper and lower bounds. Two-dimensional 
models may also be used to qualitatively evaluate the stress distribution in a ply and the variation of energy 
release rates and mixed mode ratios with delamination length. 
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TABLES 


TABLE 1. MATERIAL PROPERTIES 


IM7/8552 Unidirectional Graphite/Epoxy Prepreg 

Ell = 161.0 GPa 
V 12 = 0.32 
Gl 2 = 5.17 GPa 

E 22 = 11-38 GPa 
vi 3 = 0.32 
Gl 3 = 5.17 GPa 

£33 = 11.38 GPa 
V23 = 0.45 
G23 = 3.92 GPa 

IM7/8552 Graphite/Epoxy Plain Weave Fabric 

Ell = 71.7 GPa 
vi 2 = 0.04 
Gi 2 = 4.48 GPa 

E22 = 71.7 GPa 
vi 3 = 0.35 
Gi 3 = 4.14 GPa 

E 33 = 10.3 GPa 
V23 = 0.35 
G 23 = 4.14 GPa 

E= 1.72 GPa 

Grade 5 FM300 Adhesive 
v = 0.3 

(assumed isotropic) 
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FIGURE CAPTIONS 


Figure 1 : Skin/Flange Specimen 

(a) Specimen Configuration 

(b) Detail of typical damage patterns at the specimen comers 2 and 3 [7] 

Figure 2: Two-dimensional finite element model of skin/flange specimen 

(a) Outline of model with load and boundary conditions 

(b) Detail of model of undamaged specimen 

(c) Detail of deformed model of damaged specimen 

Figure 3: Three-dimensional generalized plane strain model of skin/flange specimen 

(a) Three-dimensional generalized plane strain model with load and boundary conditions 

(b) Detail of model of undamaged specimen 

(c) Detail of deformed model of damaged specimen 
Figure 4: Full three-dimensional model of skin/flange specimen 

(a) Three-dimensional model with load and boundary conditions 

(b) Detail of model of undamaged specimen 

(c) Detail of deformed model of damaged specimen 
Figure 5: Flange strain-load plots for quasi-static tests 
Figure 6: Typical skin strain-load plots for quasi-static tests 

Figure 7: Computed transverse tensile stress results in top 45° skin ply 

Figure 8: Computed total energy release rate for delamination growing from matrix crack in skin top 457- 
45° ply interface 

Figure 9: Computed mixed mode ratio G s /G r for delamination growing from matrix crack in skin top 
457-45° ply interface 
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(a) Specimen configuration 
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(b) Detail of typical damage patterns at specimen corners 2 and 3 [7] 
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Figure 1 . Skin/Flange Specimen 
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Static failure load: P=17.8 kN 
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(a) Outline of model with load and boundary conditions 



(b) Detail of model of undamaged specimen 



(c) Detail of deformed model of damaged specimen 
Figure 2 . Two dimensional finite element model of skin/flange specimen 


w=0 at midplane for all x and y 


v=0 at x>101.6 
P= 17.8 kN 


w=constant at parallel plane for all x and y 
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w=0 at x=0 


(a) Three-dimensional plane strain model with load and boundary conditions 
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(b) Detail of model of undamaged specimen 



(c) Detail of deformed model of damaged specimen 
Figure 3. Three-dimensional generalized plane strain model of skin/flange specimen 
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Figure 5. Flange strain-load plots for quasi-static tests. 
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Figure 6. Typical skin strain-load plot for quasi-static tests. 
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Figure 7. Computed transverse tensile stress results in top 45° skin ply. 
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Figure 8. Computed total energy release rate for delamination growing 
from matrix crack in skin top 45°/-45°ply interface 
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Figure 9. Computed mixed mode ratio G^/G ^ for delamination growing 
from matrix crack in skin top 45%45°ply interface 


